The source code for this article is available in this repository.
Note: the html version of this notebook was created with nbconvert and running the command jupyter nbconvert --to html --TagRemovePreprocessor.remove_input_tags="notvisible" notebook.ipynb. A tag "notvisible" was added to the cells that are not displayed in this rendered html
Arduino is an open-source electronics platform that provides both software and hardware to experiment with tiny machine learning, aka tinyML. An Arduino board can read inputs from sensors and turn them into output. In this project, we will use the on-board sensors for temperature, pressure and humidity. Pressure, humidity and temperature are sensors that are already integrated with Arduino. The objective is to use live readings and a model previously trained and deployed to the Arduino board to make predictions on the chance for precipitation. We will also build an iOS app and deploy it to an iPhone for acquiring the results from the board and displaying them. The app and the board do not need any Internet connection and therefore, the solution can be thought of as a companion to take to any remote region.
The code and notes in this notebook are only a summary of the work. The associated repository contains the full project.
Sub-hourly data is publicly available from the NOAA. The documentation about the data is available here. We find the Iowa State University interactive portal an easy way to understand and browse through the weather data collected at NOAA. Data can be also downloaded directly from this site. However, sub-hourly data is not provided.
The NOAA data is originally in METAR format, typically used for reporting weather information. For a list of useful acronyms and codes this format uses, the reader can look at this document, for example.
The data needed for our purpose is the timestamp, pressure, hummidity, temperature and the precipitation.
# Example of using the metar library (git+https://github.com/python-metar/python-metar.git) to convert METAR format
#Metar.Metar("63870K0J4 0J420140101125516201/01/14 12:55:31 5-MIN K0J4 011755Z AUTO 06003KT 10SM SCT005 OVC080 11/08 A3027 0 86 -400 060/03 RMK AO2 RAB04E49 SLP248 P0000 60000 T01060083 10106 20072 51006", strict=False).string().split('\n')
Metar.Metar('63870K0J4 0J420140101000012001/01/14 00:00:31 5-MIN K0J4 010500Z AUTO 03003KT 10SM -RA OVC100 07/05 A3036 -90 89 -1000 030/03 RMK AO2 RAB0459 P0000', strict=False).string().split('\n')
#Metar.Metar('63870K0J4 0J420140130094010101/30/14 09:40:31 5-MIN K0J4 301440Z AUTO 09007KT 10SM CLR M03/M06 A3034 -60 85 -2300 090/07 RMK AO2', strict=False).string().split('\n')
['station: K0J4', 'type: routine report, cycle 5 (automatic report)', 'time: Mon Apr 1 05:00:00 2024', 'temperature: 7.0 C', 'dew point: 5.0 C', 'wind: NNE at 3 knots', 'visibility: 10 miles', 'pressure: 1028.1 mb', 'weather: light rain', 'sky: overcast at 10000 feet', '1-hour precipitation: 0.00in', 'remarks:', '- Automated station (type 2)', '- RAB0459', 'METAR: 63870K0J4 0J420140101000012001/01/14 00:00:31 5-MIN K0J4 010500Z AUTO 03003KT 10SM -RA OVC100 07/05 A3036 -90 89 -1000 030/03 RMK AO2 RAB0459 P0000']
# KPHL: Philadelphia airport
# KILG: Wilmington
# KLNS: Lancaster airport
# KPNE: Philadelphia NE Airport
# KGED: Georgetown airport
years = range(2014, 2024)
months = ['{:02d}'.format(month) for month in range(1,13)]
stations = ['KILG', 'KPHL', 'KLNS', 'KPNE']
path = 'https://www.ncei.noaa.gov/data/automated-surface-observing-system-five-minute/access'
for st in stations:
for yr in years:
for mo in months:
url = f'{path}/{yr}/{mo}/asos-5min-{st}-{yr}{mo}.dat'
os.makedirs(f'../data/{st}', exist_ok=True)
urllib.request.urlretrieve(url, f'../data/{st}/{yr}{mo}.dat')
The downloaded data is preprocessed. In particular, we created the function process_file to check the validity of the data and convert it from the METAR format (see function details in the repository). The humidity is calculated from the air temperature and dew point. We use backward and forward filling of missing values with the first available non-missing value at an earlier or later timestamp.
years = range(2014, 2025)
months = ['{:02d}'.format(month) for month in range(1,13)]
stations = ['KILG', 'KPHL', 'KLNS', 'KPNE']
df = pd.DataFrame()
for st in stations:
print(f'Processing: {st}')
df_stat = pd.DataFrame()
for yr in years:
#print(f'{st} {yr}')
for mo in months:
fn = f'../data/{st}/{yr}{mo}.dat'
if os.path.exists(fn):
data = process_file(fn, )
#see https://www.ncei.noaa.gov/data/automated-surface-observing-system-five-minute/doc/asos-5min_documentation.pdf for abnormal values
data = (
data
.assign(STATION = lambda x: x['station'].astype(str))
.assign(STATION = lambda y: y['STATION'].apply(lambda x: None if x=='' or x=='None' else x)).ffill().bfill()
.assign(DATETIME = lambda x: x['datetime'].apply(lambda x: datetime.strptime(x, "%m/%d/%y %H:%M:%S")))
.sort_values('DATETIME', ascending=True)
.assign(AIR_TEMPERATURE = lambda x: pd.to_numeric(x['temperature']))
.assign(AIR_TEMPERATURE = lambda y: y['AIR_TEMPERATURE'].apply(lambda x: np.nan if x<-62 or x>54 else x).ffill().bfill()) #abnormal values
.assign(DEWPOINT = lambda x: pd.to_numeric(x['dew point']))
.assign(DEWPOINT = lambda y: y['DEWPOINT'].apply(lambda x: np.nan if x<-34 or x>30 else x).ffill().bfill()) #abnormal values
.assign(PRESSURE = lambda x: pd.to_numeric(x['pressure']))
.assign(PRESSURE = lambda y: y['PRESSURE'].apply(lambda x: np.nan if x<572.3 or x>=1066.7 else x).ffill().bfill()) #abnormal sea level pressures
.assign(RELATIVE_HUMIDITY = lambda x: calculate_relative_humidity(x['AIR_TEMPERATURE'].values, x['DEWPOINT'].values))
.assign(RELATIVE_HUMIDITY = lambda y: y['RELATIVE_HUMIDITY'].apply(lambda x: np.nan if x>100 else x).ffill().bfill()) #abnormal values
.assign(PRECIPITATION = lambda x: x['weather'].apply(lambda x: 'noprecip' if x=='' else x))
.assign(PRECIPITATION = lambda x: x['PRECIPITATION'].apply(lambda x: ('rain' in x)or('snow' in x)or('hail' in x)or('ice' in x)or('drizzle' in x)).astype(int))
.drop(['station','datetime','temperature', 'pressure', 'weather', 'dew point', 'DEWPOINT'], axis=1)
.drop_duplicates()
)
# groupby DATETIME so that we get only 1 data point for each time
data = data.groupby(['STATION','DATETIME']).agg(
{
'AIR_TEMPERATURE': 'mean',
'PRESSURE': 'mean',
'RELATIVE_HUMIDITY': 'mean',
'PRECIPITATION': 'max' # if for the same datetime we get 0 and 1 keep 1
}
).reset_index()
df_stat = pd.concat([df_stat, data], ignore_index=True)
df_stat = df_stat.sort_values(by=['STATION','DATETIME']).reset_index(drop=True)
df_stat = complete_dataframe(df_stat)
df_stat.to_parquet(f'../data/processed/{st}.pqt', index=False)
#df = pd.concat([df, df_stat], ignore_index=True)
#df = df.sort_values(by=['STATION','DATETIME']).reset_index(drop=True)
#del data
Processing: KILG Processing: KPHL Metar error Processing: KLNS Processing: KPNE
The preprocessed data is saved in Parquet format. We will add the target to predict, which we define as 1 if we observe rain in the next 30 min (6 rows, each being 5 min) else 0.
for st in stations:
df = pd.read_parquet(f'../data/processed/{st}.pqt')
display(df.isna().sum())
df = pd.read_parquet('../data/processed/KILG.pqt')
df.head()
| DATETIME | STATION | AIR_TEMPERATURE | PRESSURE | RELATIVE_HUMIDITY | PRECIPITATION | |
|---|---|---|---|---|---|---|
| 0 | 2014-01-01 00:00:31 | KILG | -3.0 | 1027.40 | 58.510 | 0 |
| 1 | 2014-01-01 00:05:00 | KILG | -3.5 | 1027.40 | 60.775 | 0 |
| 2 | 2014-01-01 00:05:31 | KILG | -4.0 | 1027.40 | 63.040 | 0 |
| 3 | 2014-01-01 00:10:00 | KILG | -4.0 | 1027.25 | 63.040 | 0 |
| 4 | 2014-01-01 00:10:31 | KILG | -4.0 | 1027.10 | 63.040 | 0 |
# add target: rain in the next 30 min (6 steps)
df['PRECIPITATION_6'] =df['PRECIPITATION'].rolling(window=6, min_periods=1).sum().shift(-6, fill_value=0).ge(1).astype(int)
# This will help us identify any data issue that needs to be addressed
df.describe().T
| count | mean | min | 25% | 50% | 75% | max | std | |
|---|---|---|---|---|---|---|---|---|
| DATETIME | 2032429 | 2018-12-18 14:36:05.884722944 | 2014-01-01 00:00:31 | 2016-06-23 12:45:00 | 2018-12-07 05:35:31 | 2021-06-22 01:50:00 | 2023-12-31 23:55:31 | NaN |
| AIR_TEMPERATURE | 2032429.0 | 13.665588 | -16.7 | 5.6 | 13.9 | 22.2 | 36.7 | 10.054358 |
| PRESSURE | 2032429.0 | 941.737224 | 580.0 | 1000.0 | 1012.9 | 1021.3 | 1045.0 | 147.319471 |
| RELATIVE_HUMIDITY | 2032429.0 | 67.701518 | 11.51 | 52.28 | 68.75 | 84.59 | 100.0 | 19.29721 |
| PRECIPITATION | 2032429.0 | 0.07857 | 0.0 | 0.0 | 0.0 | 0.0 | 1.0 | 0.269067 |
| PRECIPITATION_6 | 2032429.0 | 0.092602 | 0.0 | 0.0 | 0.0 | 0.0 | 1.0 | 0.289875 |
df.isna().sum()
DATETIME 0 STATION 0 AIR_TEMPERATURE 0 PRESSURE 0 RELATIVE_HUMIDITY 0 PRECIPITATION 0 PRECIPITATION_6 0 dtype: int64
As we can see below, the dataset is highly imbalanced with only about 9% of the rows having precipitation values equal to 1.
target_nm = 'PRECIPITATION_6'
neg, pos = np.bincount(df[target_nm])
total = neg + pos
print('Examples:\n Total: {}\n Positive: {} ({:.2f}% of total)\n'.format(
total, pos, 100 * pos / total))
class_perc = {'0': np.round(100 * neg / total, 2), '1': np.round(100 * pos / total,2)}
class_counts = {'0': neg, '1': pos}
Examples:
Total: 2032429
Positive: 188208 (9.26% of total)
fig = go.Figure()
for class_label, count in class_counts.items():
fig.add_trace(go.Bar(x=[class_label], y=[count], name=class_label))
fig.update_layout(title='Target Distribution', xaxis_title=target_nm, yaxis_title='Count')
# Add text annotations for class counts
for class_label in class_counts:
fig.add_annotation(
x=class_label, y=class_counts[class_label], text=f'{(class_perc[class_label])}%',
showarrow=True, arrowhead=4, ax=0, ay=-40
)
fig.show()
df[df[target_nm]>0].head(12)
| DATETIME | STATION | AIR_TEMPERATURE | PRESSURE | RELATIVE_HUMIDITY | PRECIPITATION | PRECIPITATION_6 | |
|---|---|---|---|---|---|---|---|
| 981 | 2014-01-02 16:55:00 | KILG | -0.8 | 1012.9 | 84.125 | 0 | 1 |
| 982 | 2014-01-02 16:55:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 983 | 2014-01-02 17:00:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 984 | 2014-01-02 17:00:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 985 | 2014-01-02 17:05:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 986 | 2014-01-02 17:05:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 987 | 2014-01-02 17:10:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 988 | 2014-01-02 17:10:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 989 | 2014-01-02 17:15:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 990 | 2014-01-02 17:15:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 991 | 2014-01-02 17:20:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 992 | 2014-01-02 17:20:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
df[981-10:981+20]
| DATETIME | STATION | AIR_TEMPERATURE | PRESSURE | RELATIVE_HUMIDITY | PRECIPITATION | PRECIPITATION_6 | |
|---|---|---|---|---|---|---|---|
| 971 | 2014-01-02 16:30:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 0 |
| 972 | 2014-01-02 16:30:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 0 |
| 973 | 2014-01-02 16:35:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 0 |
| 974 | 2014-01-02 16:35:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 0 |
| 975 | 2014-01-02 16:40:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 0 |
| 976 | 2014-01-02 16:40:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 0 |
| 977 | 2014-01-02 16:45:00 | KILG | -0.8 | 1012.9 | 84.125 | 0 | 0 |
| 978 | 2014-01-02 16:45:31 | KILG | -0.6 | 1012.9 | 81.960 | 0 | 0 |
| 979 | 2014-01-02 16:50:00 | KILG | -0.6 | 1012.9 | 81.960 | 0 | 0 |
| 980 | 2014-01-02 16:50:31 | KILG | -0.6 | 1012.9 | 81.960 | 0 | 0 |
| 981 | 2014-01-02 16:55:00 | KILG | -0.8 | 1012.9 | 84.125 | 0 | 1 |
| 982 | 2014-01-02 16:55:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 983 | 2014-01-02 17:00:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 984 | 2014-01-02 17:00:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 985 | 2014-01-02 17:05:00 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 986 | 2014-01-02 17:05:31 | KILG | -1.0 | 1012.9 | 86.290 | 0 | 1 |
| 987 | 2014-01-02 17:10:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 988 | 2014-01-02 17:10:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 989 | 2014-01-02 17:15:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 990 | 2014-01-02 17:15:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 991 | 2014-01-02 17:20:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 992 | 2014-01-02 17:20:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 993 | 2014-01-02 17:25:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 994 | 2014-01-02 17:25:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 995 | 2014-01-02 17:30:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 996 | 2014-01-02 17:30:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 997 | 2014-01-02 17:35:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 998 | 2014-01-02 17:35:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 999 | 2014-01-02 17:40:00 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
| 1000 | 2014-01-02 17:40:31 | KILG | -1.0 | 1012.9 | 86.290 | 1 | 1 |
We will use class weights to balance the dataset later during model training.
# Class weights for balancing later
class_weights = compute_class_weight('balanced', classes=np.unique(df[target_nm]), y=df[target_nm])
# Create a dictionary to pass to the fit method
class_weight = {0: class_weights[0], 1: class_weights[1]}
class_weight
{0: 0.5510264225382967, 1: 5.399422447504889}
Typically, temperature, pressure and relative humidity are good predictors for precipitation although other factors that we do not consider here can improve the final results. In fact, our objective is beyond building an extremely accurate model. In any case, we can observe below how the occurrence of precipitation can be related to a decreasing pressure and temperature, and increasing humidity.
plot_cols = ['AIR_TEMPERATURE', 'RELATIVE_HUMIDITY', 'PRESSURE', 'PRECIPITATION', 'PRECIPITATION_6']
plot_features = df[plot_cols][:1500]
plot_features.index = df['DATETIME'][:1500]
_ = plot_features.plot(subplots=True)
#plt.hist2d(df['AIR_TEMPERATURE'], df['RELATIVE_HUMIDITY'], bins=(50, 50), vmax=400)
#plt.colorbar()
#plt.xlabel('AIR_TEMPERATURE')
#plt.ylabel('RELATIVE_HUMIDITY')
#ax = plt.gca()
#ax.axis('tight')
plt.scatter(df['RELATIVE_HUMIDITY'], df['AIR_TEMPERATURE'], c=df['PRECIPITATION'], cmap='coolwarm')
plt.colorbar(label='PRECIPITATION')
plt.xlabel('RELATIVE_HUMIDITY')
plt.ylabel('AIR_TEMPERATURE')
plt.grid(True)
plt.show()
# do not shuffle
feature_nms = ['AIR_TEMPERATURE', 'RELATIVE_HUMIDITY', 'PRESSURE', target_nm]
column_indices = {name: i for i, name in enumerate(feature_nms)}
n = len(df)
train_df = df[0:int(n*0.7)][feature_nms]
val_df = df[int(n*0.7):int(n*0.9)][feature_nms]
test_df = df[int(n*0.9):][feature_nms]
# target will be integer
train_df[target_nm] = train_df[target_nm].astype(int)
test_df[target_nm] = test_df[target_nm].astype(int)
val_df[target_nm] = val_df[target_nm].astype(int)
neg, pos = np.bincount(train_df[target_nm])
total = neg + pos
print('Examples:\n Total: {}\n Positive: {} ({:.2f}% of total)\n'.format(
total, pos, 100 * pos / total))
Examples:
Total: 1422700
Positive: 135867 (9.55% of total)
neg, pos = np.bincount(test_df[target_nm])
total = neg + pos
print('Examples:\n Total: {}\n Positive: {} ({:.2f}% of total)\n'.format(
total, pos, 100 * pos / total))
Examples:
Total: 203243
Positive: 16140 (7.94% of total)
neg, pos = np.bincount(val_df[target_nm])
total = neg + pos
print('Examples:\n Total: {}\n Positive: {} ({:.2f}% of total)\n'.format(
total, pos, 100 * pos / total))
Examples:
Total: 406486
Positive: 36201 (8.91% of total)
981 1 982 1 983 1 984 1 985 1 986 1 987 1 988 1 989 1 990 1 991 1 992 1 993 1 994 1 995 1 996 1 997 1 998 1 999 1 Name: PRECIPITATION_6, dtype: int64
We preprocess the data before training the model. The input variables are standardized and we note down the means and standard deviations of the training dataset input variables for standardizing the new data coming from the Arduino sensors.
num_transformer = Pipeline(
steps=[
('scaler', StandardScaler())
]
)
preproc = ColumnTransformer(
transformers=[
('num', num_transformer, feature_nms[:-1]),
('passthr', 'passthrough', feature_nms[-1:])
]
)
ppl = Pipeline(steps=[('preproc', preproc)])
train_df = ppl.fit_transform(train_df)
test_df = ppl.transform(test_df)
val_df = ppl.transform(val_df)
cn_out = (list(preproc.named_transformers_['num'].named_steps['scaler'].get_feature_names_out()) +
feature_nms[-1:])
train_df = pd.DataFrame(train_df, columns=cn_out)
test_df = pd.DataFrame(test_df, columns=cn_out)
val_df = pd.DataFrame(val_df, columns=cn_out)
display(ppl)
Pipeline(steps=[('preproc',
ColumnTransformer(transformers=[('num',
Pipeline(steps=[('scaler',
StandardScaler())]),
['AIR_TEMPERATURE',
'RELATIVE_HUMIDITY',
'PRESSURE']),
('passthr', 'passthrough',
['PRECIPITATION_6'])]))])In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook. Pipeline(steps=[('preproc',
ColumnTransformer(transformers=[('num',
Pipeline(steps=[('scaler',
StandardScaler())]),
['AIR_TEMPERATURE',
'RELATIVE_HUMIDITY',
'PRESSURE']),
('passthr', 'passthrough',
['PRECIPITATION_6'])]))])ColumnTransformer(transformers=[('num',
Pipeline(steps=[('scaler', StandardScaler())]),
['AIR_TEMPERATURE', 'RELATIVE_HUMIDITY',
'PRESSURE']),
('passthr', 'passthrough',
['PRECIPITATION_6'])])['AIR_TEMPERATURE', 'RELATIVE_HUMIDITY', 'PRESSURE']
StandardScaler()
['PRECIPITATION_6']
passthrough
# Access the scaler from the ColumnTransformer
scaler = ppl.named_steps['preproc'].named_transformers_['num'].named_steps['scaler']
# Extract mean and standard deviation
mean_values = scaler.mean_
std_values = scaler.scale_
mean_values, std_values
(array([ 13.59036789, 68.24353611, 940.7671268 ]), array([ 10.17121142, 19.29721736, 148.35239889]))
train_df_ = train_df[feature_nms[:-1]].melt(var_name='Column', value_name='Normalized')
plt.figure(figsize=(12, 6))
ax = sns.violinplot(x='Column', y='Normalized', data=train_df_)
_ = ax.set_xticklabels(train_df.keys(), rotation=90)
We reshape the data to create a time series forecasting model in TensorFlow.
We need to define a time window where we sequence the input data pressure, humidity and temperature that is collected every 5 min for the past 30 min. The target is the precipitation label in the next 30 min, which we labeled earlier PRECIPITATION_6 (6 5-min observations). In principle, we could also use the PRECIPITATION value for each of the 5-min observations in the past 30 min as an input but we need an additional sensor connected to the Arduino to also measure it.
class WindowGenerator():
def __init__(self, input_width, label_width, shift,
train_df=train_df, val_df=val_df, test_df=test_df,
label_columns=None):
# Store the raw data.
self.train_df = train_df
self.val_df = val_df
self.test_df = test_df
# Work out the label column indices.
self.label_columns = label_columns
if label_columns is not None:
self.label_columns_indices = {name: i for i, name in
enumerate(label_columns)}
self.column_indices = {name: i for i, name in
enumerate(train_df.columns)}
# Work out the window parameters.
self.input_width = input_width
self.label_width = label_width
self.shift = shift
self.total_window_size = input_width + shift
self.input_slice = slice(0, input_width)
self.input_indices = np.arange(self.total_window_size)[self.input_slice]
self.label_start = self.total_window_size - self.label_width
self.labels_slice = slice(self.label_start, None)
self.label_indices = np.arange(self.total_window_size)[self.labels_slice]
# Given a list of consecutive inputs, the split_window method will convert them to a window of inputs
# and a window of labels.
# THIS IS MODIFIED FROM THE ORIGINAL split_window1 SO THAT IF A LABEL IS PASSED, THE INPUT WILL NOT HAVE LABEL DATA
def split_window(self, features):
inputs = features[:, self.input_slice, :]
labels = features[:, self.labels_slice, :]
if self.label_columns is not None:
labels = tf.stack(
[labels[:, :, self.column_indices[name]] for name in self.label_columns],
axis=-1)
# Exclude label columns from inputs
input_indices_to_keep = [i for i in range(inputs.shape[-1]) if i not in [self.column_indices[name] for name in self.label_columns]]
inputs = tf.gather(inputs, input_indices_to_keep, axis=-1)
# Slicing doesn't preserve static shape information, so set the shapes
# manually. This way the `tf.data.Datasets` are easier to inspect.
inputs.set_shape([None, self.input_width, None])
labels.set_shape([None, self.label_width, None])
return inputs, labels
def split_window1(self, features):
inputs = features[:, self.input_slice, :]
labels = features[:, self.labels_slice, :]
if self.label_columns is not None:
labels = tf.stack(
[labels[:, :, self.column_indices[name]] for name in self.label_columns],
axis=-1)
# Slicing doesn't preserve static shape information, so set the shapes
# manually. This way the `tf.data.Datasets` are easier to inspect.
inputs.set_shape([None, self.input_width, None])
labels.set_shape([None, self.label_width, None])
return inputs, labels
def plot(self, model=None, plot_col='', max_subplots=3, is_classification=True):
inputs, labels = self.example
plt.figure(figsize=(12, 8))
plot_col_index = self.column_indices[plot_col]
max_n = min(max_subplots, len(inputs))
for n in range(max_n):
plt.subplot(max_n, 1, n+1)
plt.ylabel(f'{plot_col} [normed]')
plt.plot(self.input_indices, inputs[n, :, plot_col_index],
label='Inputs', marker='.', zorder=-10)
if self.label_columns:
label_col_index = self.label_columns_indices.get(plot_col, None)
else:
label_col_index = plot_col_index
if label_col_index is None:
continue
plt.scatter(self.label_indices, labels[n, :, label_col_index],
edgecolors='k', label='Labels', c='#2ca02c', s=64)
if model is not None:
predictions = model(inputs)
if is_classification:
predictions = tf.cast(predictions > 0.5, dtype=tf.int32)
plt.scatter(self.label_indices, predictions[n, :, label_col_index],
marker='X', edgecolors='k', label='Predictions',
c='#ff7f0e', s=64)
if n == 0:
plt.legend()
plt.xlabel('Time [timesteps]')
# This method takes a time series DataFrame and convert it to a tf.data.Dataset
# of (input_window, label_window) pairs using the tf.keras.utils.timeseries_dataset_from_array function:
def make_dataset(self, data):
data = np.array(data, dtype=np.float32)
ds = tf.keras.utils.timeseries_dataset_from_array(
data=data,
targets=None,
sequence_length=self.total_window_size,
sequence_stride=1,
shuffle=True,
batch_size=32,)
ds = ds.map(self.split_window)
return ds
# Add properties for accessing splits as tf.data.Datasets using the make_dataset method you defined earlier.
@property
def train(self):
return self.make_dataset(self.train_df)
@property
def val(self):
return self.make_dataset(self.val_df)
@property
def test(self):
return self.make_dataset(self.test_df)
@property
def example(self):
"""Get and cache an example batch of `inputs, labels` for plotting."""
result = getattr(self, '_example', None)
if result is None:
# No example batch was found, so get one from the `.train` dataset
result = next(iter(self.train))
# And cache it for next time
self._example = result
return result
@example.setter
def example(self, new_example):
self._example = new_example
def __repr__(self):
return '\n'.join([
f'Total window size: {self.total_window_size}',
f'Input indices: {self.input_indices}',
f'Label indices: {self.label_indices}',
f'Label column name(s): {self.label_columns}'])
Here below two examples of how our WindowGenerator labels the inputs and the output coming from the dataset.
# Ex: single prediction (one output) 24 timestamps (24 5-min rows) later into the future, given 24 timestamps of history,
w1 = WindowGenerator(input_width=24,
label_width=1,
shift=24,
label_columns=feature_nms[-1:])
w1
Total window size: 48 Input indices: [ 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23] Label indices: [47] Label column name(s): ['PRECIPITATION_6']
# Ex: single prediction 1 timestamp later into the future, given 6 timestamps of history,
w2 = WindowGenerator(input_width=6,
label_width=1,
shift=1,
label_columns=feature_nms[-1:])
w2
Total window size: 7 Input indices: [0 1 2 3 4 5] Label indices: [6] Label column name(s): ['PRECIPITATION_6']
# Stack three slices, the length of the total window.
example_window = tf.stack([np.array(train_df[550:550+w2.total_window_size]),
np.array(train_df[651:651+w2.total_window_size]),
np.array(train_df[660:660+w2.total_window_size])])
example_inputs, example_labels = w2.split_window(example_window)
print('All shapes are: (batch, time, features)')
print(f'Window shape: {example_window.shape}')
print(f'Inputs shape: {example_inputs.shape}')
print(f'Labels shape: {example_labels.shape}')
All shapes are: (batch, time, features) Window shape: (3, 7, 4) Inputs shape: (3, 6, 3) Labels shape: (3, 1, 1)
w2.example = example_inputs, example_labels
w2.example
(<tf.Tensor: shape=(3, 6, 3), dtype=float64, numpy=
array([[[-1.13952679, 0.06407472, 0.55700396],
[-1.13952679, -0.06547763, 0.55700396],
[-1.13952679, -0.19502999, 0.55700396],
[-1.13952679, -0.19502999, 0.55565582],
[-1.13952679, -0.19502999, 0.55430767],
[-1.13952679, -0.19502999, 0.55430767]],
[[-1.23784349, 0.33147079, 0.53408555],
[-1.23784349, 0.33147079, 0.53408555],
[-1.23784349, 0.33147079, 0.53408555],
[-1.23784349, 0.33147079, 0.53408555],
[-1.23784349, 0.33147079, 0.53408555],
[-1.23784349, 0.33147079, 0.53408555]],
[[-1.23784349, 0.33147079, 0.53138927],
[-1.23784349, 0.33147079, 0.53138927],
[-1.23784349, 0.33147079, 0.53138927],
[-1.23784349, 0.33147079, 0.53037817],
[-1.23784349, 0.33147079, 0.52936706],
[-1.25750684, 0.41775266, 0.52936706]]])>,
<tf.Tensor: shape=(3, 1, 1), dtype=float64, numpy=
array([[[0.]],
[[0.]],
[[0.]]])>)
# this is how we can iterate of the dataset train for example
for example_inputs, example_labels in w2.train.take(1):
print(f'Inputs shape (batch, time, features): {example_inputs.shape}')
print(f'Labels shape (batch, time, features): {example_labels.shape}')
print(len(np.where(example_labels==1)[0]))
Inputs shape (batch, time, features): (32, 6, 3) Labels shape (batch, time, features): (32, 1, 1) 3
# Each element is an (inputs, label) pair.
w2.train.element_spec
(TensorSpec(shape=(None, 6, 3), dtype=tf.float32, name=None), TensorSpec(shape=(None, 1, 1), dtype=tf.float32, name=None))
We build a neural network model by searching for an optimal architecture that minimizes the Binary Focal Crossentropy loss function. This loss function seems to work well for this case with highly imbalanced labels. We also use the class weights defined earlier to mitigate the class imbalancement.
def compile_and_fit(model, window, patience=2, max_epochs=20, class_weights=None):
early_stopping = tf.keras.callbacks.EarlyStopping(monitor='val_loss', restore_best_weights=True,
patience=patience)
model.compile(#loss='binary_crossentropy',
loss = tf.keras.losses.BinaryFocalCrossentropy(gamma=2, alpha=1/class_weights[1]), # 1-alpha is weight of class 0
optimizer='adam',
metrics=['accuracy'],) # sample_weight_mode = 'temporal', ,tf.keras.metrics.F1Score()
#display(model.summary())
history = model.fit(window.train, epochs=max_epochs,
class_weight=class_weight,
validation_data=window.val,
callbacks=[early_stopping])
return history
Inputs shape (batch, time, features): (32, 1, 3) Labels shape (batch, time, features): (32, 1, 1)
We use KerasTuner to search for the optimal architecture.
## using keras tuner
tf.keras.backend.clear_session()
def model_builder(hp):
model = tf.keras.Sequential()
# Shape: (time, features) => (time*features)
model.add(tf.keras.layers.Flatten())
# tune the number of layers
hp_layers = hp.Int('layers', min_value=1, max_value=3, step=1)
for ll in range(hp_layers):
# Tune the number of units in this layer
hp_units = hp.Int(f'units_{ll}', min_value=1, max_value=32, step=1)
hp_activations = hp.Choice(f'activation_{ll}', ['relu', 'tanh', 'sigmoid'])
model.add(tf.keras.layers.Dense(units=hp_units, activation=hp_activations))
hp_dropout = hp.Float(f'dropout_{ll}', min_value=0, max_value=0.5)
model.add(tf.keras.layers.Dropout(hp_dropout))
model.add(tf.keras.layers.Dense(units=1, activation='sigmoid'))
# Add back the time dimension to match the label shape
# Shape: (outputs) => (1, outputs)
model.add(tf.keras.layers.Reshape([1, -1]))
#ThresholdingLayer(threshold=0.5) # NOTE: the label are converted to float somewhere even though they're integers
# Tune the learning rate for the optimizer
# Choose an optimal value from 0.01, 0.001, or 0.0001
hp_learning_rate = hp.Choice('learning_rate', values=[1e-3])
model.compile(optimizer=tf.keras.optimizers.Adam(learning_rate=hp_learning_rate),
loss=tf.keras.losses.BinaryFocalCrossentropy(gamma=2, alpha=1/class_weights[1]),
metrics=['accuracy'])
return model
max_epochs=100
patience = 2
# The algorithm trains a large number of models for a few epochs and carries forward
# only the top-performing half of models to the next round
tuner = kt.Hyperband(model_builder,
objective='val_loss',
max_epochs=max_epochs,
factor=3,
directory='my_dir',
project_name='rain_forecast')
early_stopping = tf.keras.callbacks.EarlyStopping(monitor='val_loss',
patience=patience)
tuner.search(wide_window.train, validation_data=wide_window.val,
callbacks=[early_stopping], class_weight=class_weight)
Trial 254 Complete [00h 04m 18s] val_loss: 0.1372850090265274 Best val_loss So Far: 0.11420851200819016 Total elapsed time: 17h 18m 15s
tuner.results_summary(num_trials=1)
Results summary Results in my_dir/rain_forecast Showing 1 best trials Objective(name="val_loss", direction="min") Trial 0142 summary Hyperparameters: layers: 3 units_0: 18 activation_0: relu dropout_0: 0.03492223941614769 learning_rate: 0.001 units_1: 8 activation_1: relu dropout_1: 0.31268099900849095 units_2: 4 activation_2: relu dropout_2: 0.2108992311717759 tuner/epochs: 34 tuner/initial_epoch: 12 tuner/bracket: 4 tuner/round: 3 tuner/trial_id: 0131 Score: 0.11420851200819016
Once we find the optimal neural network architecture, we train our model, we measure its performance and save it.
tf.keras.backend.clear_session()
dense = tf.keras.Sequential([
# Shape: (time, features) => (time*features)
tf.keras.layers.Flatten(),
tf.keras.layers.Dense(units=18, activation='relu'),
tf.keras.layers.Dropout(0.03),
tf.keras.layers.Dense(units=8, activation='relu'),
tf.keras.layers.Dropout(0.3),
tf.keras.layers.Dense(units=4, activation='relu'),
tf.keras.layers.Dropout(0.2),
#tf.keras.layers.BatchNormalization(),
#tf.keras.layers.Dense(units=3, activation='sigmoid'),
tf.keras.layers.Dense(units=1, activation='sigmoid'),
# Add back the time dimension to match the label shape
# Shape: (outputs) => (1, outputs)
tf.keras.layers.Reshape([1, -1]),
])
print('Input shape:', wide_window.example[0].shape)
print('Output shape:', dense(wide_window.example[0]).shape)
history = compile_and_fit(dense, wide_window, max_epochs=100, patience=10, class_weights=class_weight)
Input shape: (32, 1, 3) Output shape: (32, 1, 1) Epoch 1/100 44460/44460 [==============================] - 44s 978us/step - loss: 0.1391 - accuracy: 0.6802 - val_loss: 0.1219 - val_accuracy: 0.7672 Epoch 2/100 44460/44460 [==============================] - 47s 1ms/step - loss: 0.1376 - accuracy: 0.6894 - val_loss: 0.1232 - val_accuracy: 0.7647 Epoch 3/100 44460/44460 [==============================] - 48s 1ms/step - loss: 0.1373 - accuracy: 0.6908 - val_loss: 0.1190 - val_accuracy: 0.7651 Epoch 4/100 44460/44460 [==============================] - 51s 1ms/step - loss: 0.1370 - accuracy: 0.6927 - val_loss: 0.1171 - val_accuracy: 0.7960 Epoch 5/100 44460/44460 [==============================] - 54s 1ms/step - loss: 0.1369 - accuracy: 0.6960 - val_loss: 0.1249 - val_accuracy: 0.7613 Epoch 6/100 44460/44460 [==============================] - 50s 1ms/step - loss: 0.1370 - accuracy: 0.6960 - val_loss: 0.1221 - val_accuracy: 0.7834 Epoch 7/100 44460/44460 [==============================] - 50s 1ms/step - loss: 0.1368 - accuracy: 0.6977 - val_loss: 0.1232 - val_accuracy: 0.7679 Epoch 8/100 44460/44460 [==============================] - 104s 2ms/step - loss: 0.1369 - accuracy: 0.6946 - val_loss: 0.1239 - val_accuracy: 0.7810 Epoch 9/100 44460/44460 [==============================] - 45s 1ms/step - loss: 0.1367 - accuracy: 0.6950 - val_loss: 0.1248 - val_accuracy: 0.7789 Epoch 10/100 44460/44460 [==============================] - 51s 1ms/step - loss: 0.1367 - accuracy: 0.6947 - val_loss: 0.1195 - val_accuracy: 0.7852 Epoch 11/100 44460/44460 [==============================] - 53s 1ms/step - loss: 0.1369 - accuracy: 0.6964 - val_loss: 0.1229 - val_accuracy: 0.7664 Epoch 12/100 44460/44460 [==============================] - 52s 1ms/step - loss: 0.1367 - accuracy: 0.6967 - val_loss: 0.1225 - val_accuracy: 0.7680 Epoch 13/100 44460/44460 [==============================] - 53s 1ms/step - loss: 0.1367 - accuracy: 0.6972 - val_loss: 0.1266 - val_accuracy: 0.7793 Epoch 14/100 44460/44460 [==============================] - 54s 1ms/step - loss: 0.1368 - accuracy: 0.6985 - val_loss: 0.1264 - val_accuracy: 0.7822
dense.summary()
Model: "sequential"
_________________________________________________________________
Layer (type) Output Shape Param #
=================================================================
flatten (Flatten) (None, 3) 0
dense (Dense) (None, 18) 72
dropout (Dropout) (None, 18) 0
dense_1 (Dense) (None, 8) 152
dropout_1 (Dropout) (None, 8) 0
dense_2 (Dense) (None, 4) 36
dropout_2 (Dropout) (None, 4) 0
dense_3 (Dense) (None, 1) 5
reshape (Reshape) (None, 1, 1) 0
=================================================================
Total params: 265 (1.04 KB)
Trainable params: 265 (1.04 KB)
Non-trainable params: 0 (0.00 Byte)
_________________________________________________________________
From the performance graphs below, it seems that our model underfits the data. However, we should keep in mind that we don't have many data exemplars with positive labels and during the splitting training/ validation/ test samples positive labels exemplars were randomly assigned to the splits. Some precipitation occurrences might have "unique" associated temperature, humidity, pressure values that are only see by one of the splits.
loss_train = history.history['loss']
loss_val = history.history['val_loss']
acc_train = history.history['accuracy']
acc_val = history.history['val_accuracy']
epochs = range(1, len(history.history['loss']) + 1)
def plot_train_val_history(x, y_train, y_val, type_txt):
plt.figure(figsize = (10,7))
plt.plot(x, y_train, 'g', label='Training'+type_txt)
plt.plot(x, y_val, 'b', label='Validation'+type_txt)
plt.title('Training and Validation'+type_txt)
plt.xlabel('Epochs')
plt.ylabel(type_txt)
plt.legend()
plt.show()
plot_train_val_history(epochs,
loss_train, loss_val,
"Loss")
plot_train_val_history(epochs,
acc_train, acc_val,
"Accuracy")
# when getting batches during training, they're shuffled.
# Make sure here to retrieve the same order for metrics evaluation
# This needs to be fixed for multi-outputs
def get_labels(window, model):
train_inputs = []
train_labels = []
test_inputs = []
test_labels = []
val_inputs = []
val_labels = []
for i, j in window.train:
train_inputs.append(i)
train_labels.append(j)
train_inputs = tf.concat(train_inputs, axis=0)
train_labels = tf.concat(train_labels, axis=0).numpy().ravel()
y_trainpreds = (model.predict(train_inputs) > 0.5).astype("int32").ravel()
for i, j in window.val:
val_inputs.append(i)
val_labels.append(j)
val_inputs = tf.concat(val_inputs, axis=0)
val_labels = tf.concat(val_labels, axis=0).numpy().ravel()
y_valpreds = (model.predict(val_inputs) > 0.5).astype("int32").ravel()
for i, j in window.test:
test_inputs.append(i)
test_labels.append(j)
test_inputs = tf.concat(test_inputs, axis=0)
test_labels = tf.concat(test_labels, axis=0).numpy().ravel()
y_testpreds = (model.predict(test_inputs) > 0.5).astype("int32").ravel()
return y_trainpreds, train_labels, y_valpreds, val_labels, y_testpreds, test_labels
y_trainpreds, train_labels, y_valpreds, val_labels, y_testpreds, test_labels = get_labels(wide_window, dense)
44460/44460 [==============================] - 13s 287us/step 12703/12703 [==============================] - 4s 278us/step 6352/6352 [==============================] - 2s 275us/step
# training set performance
model_performance(y_trainpreds, train_labels, plot=False)
Accuracy: 0.78 Recall: 0.268 Precision: 0.752 F-score: 0.395
# test set performance
model_performance(y_testpreds, test_labels, plot=False)
Accuracy: 0.788 Recall: 0.225 Precision: 0.683 F-score: 0.338
We can finally look at the ROC, which is not great.
<matplotlib.legend.Legend at 0x30c5c9250>
dense.save("./model/rain_forecast")
TensorFlow Lite for microcontrollers is designed to run machine learning models on microcontrollers and other devices with only a few kilobytes of memory. A "hello world" example for Arduino is available here.
Typically, post-training integer 8-bit quantization is the main ingredient that make the trained model suitable for inference on devices with reduced memory computational capabilities such as an Arduino board. This widely adopted technique applies the quantization after training and converts the 32-bit floating-point weights into 8-bit integer values. However, we do not need to quantize the model we built due to its already small size. Quantizing the model will significantly reduce its predictive performance.
TF_MODEL = "./model/rain_forecast"
TFL_MODEL_FILE = "./model_tflite/rain_model.tflite"
is_quantize = False
# use some test data for calibration
# generating a representative dataset is essential as it helps reduce the risk of an accuracy drop
# when the model is converted to an 8-bit format. In fact, the converter uses a set of samples to
# find out the distributions of the floating-point values and estimate the optimal quantization parameters.
def representative_dataset(num_sample=500):
data = wide_window.test
for i_value, _ in data.batch(1).take(100):
i_value_f32 = tf.dtypes.cast(i_value, tf.float32)
yield [i_value_f32]
converter = tf.lite.TFLiteConverter.from_saved_model(TF_MODEL)
#converter = tf.lite.TFLiteConverter.from_keras_model(dense,)
if is_quantize:
converter.optimizations = [tf.lite.Optimize.DEFAULT]
converter.target_spec.supported_ops = [tf.lite.OpsSet.TFLITE_BUILTINS_INT8]
converter.inference_input_type = tf.int8
converter.inference_output_type = tf.int8
#converter.inference_input_type = tf.float32
#converter.inference_output_type = tf.float32
converter.representative_dataset = representative_dataset
# convert to TFlite format
tflite_model = converter.convert()
# save TFLite model
open(TFL_MODEL_FILE, "wb").write(tflite_model)
2024-04-08 16:53:16.017472: W tensorflow/compiler/mlir/lite/python/tf_tfl_flatbuffer_helpers.cc:378] Ignored output_format.
2024-04-08 16:53:16.017647: W tensorflow/compiler/mlir/lite/python/tf_tfl_flatbuffer_helpers.cc:381] Ignored drop_control_dependency.
2024-04-08 16:53:16.017940: I tensorflow/cc/saved_model/reader.cc:83] Reading SavedModel from: ./model/rain_forecast
2024-04-08 16:53:16.019138: I tensorflow/cc/saved_model/reader.cc:51] Reading meta graph with tags { serve }
2024-04-08 16:53:16.019144: I tensorflow/cc/saved_model/reader.cc:146] Reading SavedModel debug info (if present) from: ./model/rain_forecast
2024-04-08 16:53:16.021035: I tensorflow/compiler/mlir/mlir_graph_optimization_pass.cc:388] MLIR V1 optimization pass is not enabled
2024-04-08 16:53:16.022268: I tensorflow/cc/saved_model/loader.cc:233] Restoring SavedModel bundle.
2024-04-08 16:53:16.060164: I tensorflow/cc/saved_model/loader.cc:217] Running initialization op on SavedModel bundle at path: ./model/rain_forecast
2024-04-08 16:53:16.072131: I tensorflow/cc/saved_model/loader.cc:316] SavedModel load for tags { serve }; Status: success: OK. Took 54193 microseconds.
2024-04-08 16:53:16.086550: I tensorflow/compiler/mlir/tensorflow/utils/dump_mlir_util.cc:269] disabling MLIR crash reproducer, set env var `MLIR_CRASH_REPRODUCER_DIRECTORY` to enable.
Summary on the non-converted ops:
---------------------------------
* Accepted dialects: tfl, builtin, func
* Non-Converted Ops: 13, Total Ops 26, % non-converted = 50.00 %
* 13 ARITH ops
- arith.constant: 13 occurrences (f32: 8, i32: 5)
(f32: 4)
(f32: 1)
(i32: 1)
(f32: 2)
(i32: 1)
(i32: 1)
4896
size_tfl_model = len(tflite_model)
print(size_tfl_model, "bytes")
4896 bytes
Since we did not quantize the model, the TFLite version should retain the original model's performance.
from tqdm import tqdm
# TFLite interpreter
tflite_interpreter = tf.lite.Interpreter(
model_path=TFL_MODEL_FILE,
experimental_op_resolver_type=tf.lite.experimental.OpResolverType.BUILTIN_REF,
)
tflite_interpreter.allocate_tensors()
input_details = tflite_interpreter.get_input_details()[0]
output_details = tflite_interpreter.get_output_details()[0]
input_shape = np.array(input_details.get('shape'))
print(f'Expected input shape: {input_shape}')
input_index = input_details['index']
output_index = output_details['index']
# prepare data to test
train_inputs = []
train_labels = []
for i, j in wide_window.train:
train_inputs.append(i)
train_labels.append(j)
train_inputs = tf.concat(train_inputs, axis=0)
train_labels = tf.concat(train_labels, axis=0).numpy().ravel()
print(f'Train input shape: {train_inputs.shape}')
print(f'Train output shape: {train_labels.shape}')
##
# generate predictions
train_predictions = np.empty(tf.shape(train_inputs)[0], dtype=np.float32)
print(f'Expected prediction shape: {train_predictions.shape}')
def tflite_predict(input_shape, interpreter, input_value, input_index, output_index):
input_data = np.reshape(input_value, input_shape)
tflite_interpreter.set_tensor(input_index, input_data)
interpreter.invoke()
tflite_output = tflite_interpreter.get_tensor(output_index)
return np.reshape(tflite_output, -1)[0]
for i in tqdm(range(train_inputs.shape[0])):
train_predictions[i] = tflite_predict(
input_shape,
tflite_interpreter,
train_inputs[i],
input_index,
output_index,
)
train_predictions = (train_predictions>0.5).astype(int)
Expected input shape: [1 1 3] Train input shape: (1422699, 1, 3) Train output shape: (1422699,) Expected prediction shape: (1422699,)
100%|██████████| 1422699/1422699 [01:58<00:00, 11986.28it/s]
As expected, the performance remain the same after converting the TensorFlow model to a TFLite model.
model_performance(train_labels, train_predictions)
Accuracy: 0.788 Recall: 0.744 Precision: 0.275 F-score: 0.401
The model binary data need to be converted to a constant C-byte array. The steps we need to follow are in the comments below. Since the model will be consumed by the Arduino processor, we will save the model directly into the folder that we set up later in the Arduino IDE for deployment.
# The xxd tool is a hex dump utility that can be used to create a hexadecimal representation of binary data. xxd is a native command on Unix
# see https://github.com/tensorflow/tflite-micro-arduino-examples/blob/main/examples/hello_world/hello_world.ino (cpp and h for model requirements)
!xxd -i "./model_tflite/rain_model.tflite" > ./model_tflite/model.cpp
# NOTE: !sed -i 's/unsigned char/const unsigned char/g' model.h (original sybtax)
# the syntax below is sometimes required on macOS to indicate that an in-place edit is being performed,
# and it might resolve the "invalid command code" error.
#!sed -i '' 's/unsigned char/const unsigned char/g' ./model_tflite/model.cpp
#!sed -i '' 's/const/alignas(8) const/g' ./model_tflite/model.h
!cp ./model_tflite/model.cpp ../arduino/iosWeather/
# then add/ edit:
# #include "model.h"
# alignas(16) const unsigned char __model_tflite_rain_model_tflite[] = ...
# const int __model_tflite_rain_model_tflite_len = sizeof(__model_tflite_rain_model_tflite);
The deployment to the Arduino board is performed through the Arduino IDE that needs to be installed. We are not covering the development of the code needed to read the sensor data. The reader can refer to the repository for additional information and code (folder: arduino/iosWeather).

The Arduino board we use in this project is the Nano 33 BLE. Pressure, humidity and temperature are sensors already integrated with this Arduino. For information about the sensors, the user can look here.

Once we open the Arduino IDE, we create a new sketch (we call it iosWeather) and add the required libraries to it (Sketch > Include Librariy > Manage Libraries for a MacOS). In particular, we add the libraries to read from the sensors HS3003 (humidity, temperature) and LPS22HB (pressure). This can be done by searching the libraries in the Library Manager (see earlier snapshot).
We also add the model (model.h and model.cpp) to the sketch (Sketch > Add File).
Note: you need to install the TFLite libraries for Arduino (see the tutorials here and the examples). The TFlite libraries are manually installed into the Arduino library folder. The location for this folder varies by operating system, but typically it is in ~/Arduino/libraries on Linux, ~/Documents/Arduino/libraries/ on MacOS, and My Documents\Arduino\Libraries on Windows.
git clone https://github.com/tensorflow/tflite-micro-arduino-examples Arduino_TensorFlowLite
cd Arduino_TensorFlowLite
git pull
We can now use TFLite in our Arduino sketch. All the tensors the model consumes must be allocated, such as the input, output, and intermediate ones. The reason is that the TensorFlow Lite Interpreter does not allocate any memory at runtime.
The main .ino file must be written in C language.

and has two main functions

a setup function to initialize sensors, Bluetooth, arrays, tensors, and a loop function that streams the sensor data. We let the reader study the files in the repository to have more details about the structure.
With the Arduino board connected via USB to the computer, we can compile and upload the sketch. Note: uploading a new sketch or an empty sketch will replace the current sketch.
ADD IMAGE
Once the model and the Arduino script are deployed to the board, we can start testing them to validate the Arduino results match the results we obtain from running the model in our notebook.
For example, we take some rondom batches and shape the associated input and output arrays so that they match the arrays prepared in the Arduino script. Also, we run both the TensorFlow and TFLite models and save the predictions.
After we replace the actual sensor reading with the "simulated" readings we copy/ paste from here into the Arduino script, the predictions of the model deployed to the board should return the same results.
# Find examples with label==1
for i_value, o_value in wide_window.test.batch(1).take(1):
print(i_value)
print(o_value)
o_value.numpy()[0,0,:,:]
i_value.numpy()[0,0,:,:]
array([[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.12958941, 0.37843433],
[ 1.6202374 , -0.08279514, 0.37843433],
[ 1.6202374 , -0.12958941, 0.37843433],
[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.23084913, 0.37843433],
[ 1.6202374 , -0.2853146 , 0.37843433],
[ 1.644923 , -0.2732964 , 0.37843433],
[ 1.6696086 , -0.26127818, 0.37843433],
[ 1.6696086 , -0.26127818, 0.37843433]], dtype=float32)
import numpy as np
# Your input data
test_arr = [[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.12958941, 0.37843433],
[ 1.6202374 , -0.08279514, 0.37843433],
[ 1.6202374 , -0.12958941, 0.37843433],
[ 1.6202374 , -0.17638367, 0.37843433],
[ 1.6202374 , -0.23084913, 0.37843433],
[ 1.6202374 , -0.2853146 , 0.37843433],
[ 1.644923 , -0.2732964 , 0.37843433],
[ 1.6696086 , -0.26127818, 0.37843433],
[ 1.6696086 , -0.26127818, 0.37843433]]
test_arr = np.array(test_arr).reshape(1,12,3)
test_arr = tf.cast(test_arr, tf.float32)
dense.predict(test_arr)
1/1 [==============================] - 0s 21ms/step
array([[[0.25992578]]], dtype=float32)
tflite_predict(
input_shape,
tflite_interpreter,
test_arr,
input_index,
output_index,
)
0.2599258
# Generate code to copy and paste into the Arduino script
# We assume the sensors measure the values we randomly selected from the batches above.
code_lines = []
for i, sublist in enumerate(test_arr):
for j, subsublist in enumerate(sublist):
for k, value in enumerate(subsublist):
code_line = f'inputs_arr[{i}][{j}][{k}] = {value};'
code_lines.append(code_line)
# Print generated code
for line in code_lines:
print(line)
inputs_arr[0][0][0] = 1.6202373504638672; inputs_arr[0][0][1] = -0.17638367414474487; inputs_arr[0][0][2] = 0.37843433022499084; inputs_arr[0][1][0] = 1.6202373504638672; inputs_arr[0][1][1] = -0.17638367414474487; inputs_arr[0][1][2] = 0.37843433022499084; inputs_arr[0][2][0] = 1.6202373504638672; inputs_arr[0][2][1] = -0.17638367414474487; inputs_arr[0][2][2] = 0.37843433022499084; inputs_arr[0][3][0] = 1.6202373504638672; inputs_arr[0][3][1] = -0.12958940863609314; inputs_arr[0][3][2] = 0.37843433022499084; inputs_arr[0][4][0] = 1.6202373504638672; inputs_arr[0][4][1] = -0.0827951431274414; inputs_arr[0][4][2] = 0.37843433022499084; inputs_arr[0][5][0] = 1.6202373504638672; inputs_arr[0][5][1] = -0.12958940863609314; inputs_arr[0][5][2] = 0.37843433022499084; inputs_arr[0][6][0] = 1.6202373504638672; inputs_arr[0][6][1] = -0.17638367414474487; inputs_arr[0][6][2] = 0.37843433022499084; inputs_arr[0][7][0] = 1.6202373504638672; inputs_arr[0][7][1] = -0.23084913194179535; inputs_arr[0][7][2] = 0.37843433022499084; inputs_arr[0][8][0] = 1.6202373504638672; inputs_arr[0][8][1] = -0.2853145897388458; inputs_arr[0][8][2] = 0.37843433022499084; inputs_arr[0][9][0] = 1.6449229717254639; inputs_arr[0][9][1] = -0.27329638600349426; inputs_arr[0][9][2] = 0.37843433022499084; inputs_arr[0][10][0] = 1.6696085929870605; inputs_arr[0][10][1] = -0.2612781822681427; inputs_arr[0][10][2] = 0.37843433022499084; inputs_arr[0][11][0] = 1.6696085929870605; inputs_arr[0][11][1] = -0.2612781822681427; inputs_arr[0][11][2] = 0.37843433022499084;
Here we generate the Universally Unique Identifier (UUID) numbers that identify the service characteristics provided by the Bluetooth device. In the Arduino file, we can see the assignment of these UUIDs to the different characteristics we want to read and notify via the Bluetooth service.
BLEService weatherService("3c2312f7-7c24-4104-85ad-85a5a039d3cd");
BLEByteCharacteristic weatherCharacteristic( "4c40dfe1-7044-4058-8d53-24372c4e2e08", BLERead | BLENotify );
BLEFloatCharacteristic tempCharacteristic( "aaa691a2-3713-46ed-8d85-3a120d760175", BLERead | BLENotify );
BLEFloatCharacteristic humCharacteristic( "d85e531e-0d16-4e39-902c-417b100e867e", BLERead | BLENotify );
Similarly, the iosApp will have the same UUIDs used for the connectivity to the Arduino board via Bluetooth.
import uuid
uuid.uuid4()
UUID('d85e531e-0d16-4e39-902c-417b100e867e')
The app is developed using Xcode, a comprehensive suite of tools and resources that allow for designing, coding, and debugging iOS applications. When you open Xcode for the first time, it will ask what platform we want to develop for and download the required libraries.

We can now select “Create a new project” from the welcome screen and choose the app template from the iOS section. The reader can find the project in the repository (iosapp/ioasWeather). The code is in Swift language. SwiftUI is also available to help design the user interface.

!jupyter nbconvert --to html --TagRemovePreprocessor.remove_input_tags="notvisible" tflite-on-arduino.ipynb --output ../../dsprojects/posts/04_15_2024_tflite_on_arduino.html
[NbConvertApp] Converting notebook tflite-on-arduino.ipynb to html [NbConvertApp] Writing 4756648 bytes to ../../dsprojects/posts/04_15_2024_tflite_on_arduino.html
Copyright (c) [2024] [Alessio Tamburro]
Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to deal in the Software without restriction, including without limitation the rights to use, copy, modify, merge, publish, distribute, sublicense, and/or sell copies of the Software, and to permit persons to whom the Software is furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.